Time auto-correlation function and Green-kubo formula: A study on disordered harmonic chain 
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We have considered heat conduction in a one-dimensional mass disordered harmonic chain of N particles 
connected to two Langevin type reservoirs at different temperatures. An exact expression for the boundary heat 
current-current auto-correlation function in the non-equilibrium steady state (NESS) is obtained in terms of 
non-equilibrium phonon Green's functions. The time integral of the correlation function gives expected result, 
both in non-equilibrium as well as equilibrium cases. Using the form of this correlation function we show that 
asymptotic system size dependence of current fluctuation in NESS for a mass disordered harmonic chain is N~ a 
for different boundary conditions. For free and fixed boundary conditions we get a = 1/2 and 3/2 respectively, 
while for pinned case the fluctuation decays exponentially with system size. 

PACS numbers: 



I. Introduction 

Time correlation functions are useful quantities in the study 
of transport processes. They are related to various transport 
coefficients. For example, the diffusion constant of a Brown- 
ian particle is given by the integral of the equilibrium velocity- 
velocity time auto-correlation function. Similarly the friction 
coefficient of an over-damped particle is also related to the 
time correlation function of the instantaneous force experi- 
enced by the particle. Let us consider a stochastic process de- 
scribed by the vector x(i). Then the time correlation function 
of any quantity A(t) = A(x(t)) is defined as: (A(t)A(t')) 
where (} represents the average over initial conditions and 
trajectories. In terms of phase space variables, (A(t)A(t')) 
is given by 

(A(t)A(t>)} = JdxJ dyA(x)A(y)W(x,t;y,t')P(y,t')(l) 

where, P(y, t') is the probability of y at time t' and W is the 
transition probability from y to x in time t — t'. In general the 
equilibrium time correlation function of some quantity is re- 
lated to the response of a system to small perturbations. These 
relations are called Green-Kubo formula (GK) [1,2]. 

For the case of heat transport the GK formula relates the re- 
sponse of a system to a small temperature gradient to the equi- 
librium heat current auto correlation function. The response 
to temperature gradient defines the thermal conductivity k and 
the GK formula gives: 



Km lim = — -j 



dt(J(t)J(0)) 



(2) 



where J(t) is the heat current through the system at time t 
and L is the linear dimension of a d-dimensional system. In 
Eq. (2) the order of the limits is very important. Although this 
is a very useful formula, there are some difficulties associated 
with this formula. The formula in Eq. (2) is not applicable to 
small mesoscopic structures. Also in case of anomalous trans- 
port, which occurs in many low dimensional systems, the ther- 
mal conductivity diverges [3, 4]. In such cases it is not possi- 
ble to take the limits as in Eq. (2). There are various deriva- 
tions of this formula [5, 6]. Recently we have derived a for- 
mula similar to Eq. (2) for open systems, which is applicable 



to systems of arbitrary size in any dimensions [7]. This deriva- 
tion uses Fokker-planck description of stochastic systems and 
hence is only applicable for those currents, which can be ex- 
pressed in terms of phase space variables (e.g. currents inside 
the bulk of the system). Since boundary currents naturally in- 
volve noises explicitly, derivation given in [7] is not be appli- 
cable for them. General expectation is, for boundary currents 
also one can proof a open finite system GK formula as given 
in [7]. In this paper we explicitly calculate boundary current- 
current auto correlation function in the context of heat trans- 
port for a finite mass disordered harmonic chain in NESS and 
show that integration of the equilibrium correlation function 
gives the NESS current. 

There are few examples where exact time auto-correlation 
functions in equilibrium state have been obtained for many- 
particle systems. For Hamiltonian systems some examples 
of exact calculations are velocity auto-correlation function for 
ordered harmonic lattices [8] and for a one dimensional gas 
of elastically colliding hard rods [9]. Recently authors of 
[10] have shown explicitly that integration of the heat current 
auto-correlation function gives the current in non-equilibrium 
steady state for a two particle harmonic system. In this paper 
we obtain an exact expression for the time auto-correlation 
function for boundary heat current in the NESS for mass disor- 
dered harmonic chains of arbitrary length, expressed in terms 
of the non-equilibrium Green's functions. We show that it sat- 
isfies the GK formula derived in [7]. Using this correlation 
function we also calculate the asymptotic system size scaling 
of fluctuations in current in NESS. 

The paper is organised as follows. In sec. (II) we give the 
description of the model, define some relevant quantities and 
notations and calculate the current in the NESS. In sec. (Ill) 
we present the calculation of the time correlation function. 
In sec. (IV) we discuss our results and finally in sec.(V) we 
conclude. 



II. Definition of model 

We consider a chain of oscillators of N particles described 
by the Hamiltonian H : 
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N 1 1 

i=i 



N-l 
1=1 



k(xi + i - xiY 



+ 2^'( x \ + x n) i 



(3) 



where x t are displacements of the particles about their equilib- 
rium positions, k, ko are the inter-particle and on-site spring 
constants respectively, and to/ is mass of the I th particle. 
k! is the spring constant of the potentials at the boundaries. 
For different values of k! and ko we get different boundary- 
conditions (BCs). If k' and ko both are zero we get free BC, 
otherwise we get fixed BC (k' ^ and ko = 0) and pinned 
case (ko ^ 0). The particles 1 and N are connected to two 
white noise heat baths of temperatures Tl and T R respec- 
tively. The equation of motion of the I th particle is given by 
[11] 

miii = —k(2xi - xi-! - xi+i) - k Q xi 
—Si,i[(k' - k)xi + 7l±i - til] 
Si,N[{k' - k)xi + j R x N - r] R } 
where I — 1,2...N and x = = 0(4) 

where r]L, R (t) are Gaussian noise terms with zero mean and 
related to the dissipative terms with these relations 

(VL,R{t)VL,R(t')) = 2lL,RT L , R 5{t - t') 

(VL(t)m(t')) = 0, <rjL,fl(t)> = (5) 

(In this paper we have set Kb = 1.) To define the local energy 
current inside the chain we first define the local energy density 
associated with the / th particle (or energy at the lattice site I) 
as follows: 



2to-i 

ei = ik l + ^ L + \ [i * Xi - i ~ Xi)2 + {xi ~ Xi+i)2 



for I = 2,3..JV- 1 



Pn 
2m n 



N 



'N 



k 



., + - x n) 2 ■ (6) 



Using this energy density we write a continuity equation, from 
which we get two instantaneous currents ]l and j R which are 
flowing from the left and right reservoirs into the system re- 
spectively. These currents are given by [3, 4] 



Mt) = -7£*i(*)+»te(t)£i(t) 



and j R (t) 



-lRX 2 N {t) + rj R (t)x N {t) 



(7) 



In order to obtain the steady state properties we have to find 
out the steady state solution of the Eq. (4). For that we write 
Eq. (4) in Matrix form as: 

MX + TX + $X = r](t), (8) 

where, X, r\ are column vectors with elements [X] T = 

(xi,x 2 , ■■■■xn), [n} T = {Vl,0, ■•••0,T?fl) and T is a iVx AT ma- 
trix with only non-vanishing elements [T] 11 = jl, [F] nn = 



1r- W\nx.n represents a tridiagonal matrix with elements 
[12] 

= {k + k' + k )5^ m - kSl.m-! for I = 1 

= -k6i, m -i + (2k + k )5 Lm - k5^ m+ i 

for 2 < I < N - 1 
= (k + k' + k )5 u - k5 Lm+l for I = N , (9) 

and Mi m = miSi m where to; is chosen uniformly from the 
range [1 — A, 1 + A]. If Q + (t) denotes the Green's function 
of the entire system then Q + (t) satisfies 



Mg+(t) + vg+{t) + <s>g+(t) = s(t)i , 



(10) 



It is easy to verify that Q + (t) = G(t)Q(t) where G(t) satisfies 
the homogeneous equation 



MG + TG + $G = , 



(11) 



with the initial conditions G(0) = 0, G(0) = M -1 . Here 
Q(t) is the Heaviside function. Assuming that the heat baths 
have been switched on at t = — oo we write the steady state 
solution of Eq. (8) as 

X(t) = I dt'G(t - t')r](t'). (12) 

J — oo 

For equilibration we require that G(t) — > as t — > oo. From 
Eq.(12), we get 



Xl (t) = f 

J —I 



Guit-h^Lih) 



+G 1N (t-t 1 )r] R (t 1 ) 



(13) 



Next we calculate (j L ) in the NESS. Here (...) denotes the 
average over the noise variables and rj R (t). From now 

we denote IJl) by j. Putting x\(t) from Eq. (13) in the ex- 
pression of ]l (t) in Eq. (7) and using the noise correlation in 
Eq. (5) we get : 



3 



-1L 



f dti f 

J —oo J —i 



dU 



G n (t - ti)Gn(t - t 2 ) 



*(VL(tl)VL(t2)) + G 1N (t - h)G 1N (t - t 2 ) 



x{VR(h)VR(t2)) 



+ 



f 

J — ( 



dhGuit - h)(r] L (t)r] L {h)) 



Gn(0) - ( 7L T L Ai(0) + ~f R T R A N (0)) 



(14) 



where we have used the definition 



Ai(t)= dt'Gu(t + t')G u (t') V t. (15) 
Jo 

We now note the following identity (for proof see Appendix 
A) 



Tl-Ai (t) + ^y R Apf{i) 



Gn(t) 



(16) 



which can be obtained from Eqs.(15,ll). Using this in 
Eq. (14) we get 

3 =2jlJr(T l -T r )A n (0) . (17) 

If we go to the frequency u> space using the following defini- 
tion 



g+h = / dtG(ty 

Jo 



we can identify that 



and 



G+ (uj) = - Muj' 2 + iujT + $ 



(18) 



(19) 



(20) 



With this identification we see that the expression given in 
Eq. (17) reduces to the form 



(Tl - T R ) 
2tt 



duj T(tu) 



where 



T(uj) = 4 7L 7_r oj 2 \Gf N (w)\ 



(21) 



(22) 



is the transmission coefficient for frequency uj. The above 
expression for the current j is seen to be identical to the well- 
known expression for the current given in [13, 14]. 
In the next section we proceed to obtain the time auto- 
correlation function Cat^, t') defined as: 



Gat(M') = (3L(t)jL(t')} ~ (j L f 



(23) 



in the NESS. The subscript AT represents the difference be- 
tween the temperature at the two ends i.e. AT = Tl — T R . In 
the stationary state (ii(t)ii(t')) will t> e a function of \t — t'\ 
only. Hence we set t' — 0. If we take AT = in the expres- 
sion of Gat(^) we get the equilibrium auto-correlation which 
is denoted by Co(t) and we show that integral of Co(t) is re- 
lated to the average current (j'l), whereas integral of Gat(^) 
is related to its fluctuations in the NESS. 



III. Calculation of auto-correlation function 

Using the forms of jl from Eq. (7) we write current current 
auto-correlation (jL(t)jL(0)) as: 

0'i(*)ii(0)) = Jli + J L 2 + J L 2 + JL4 , 
where 

Jli = 7i<i?(*)i?(0)>, 

J L 2 = -7L<»?z(t)ii(t)±i(0)}, 

J L3 = -7L<T?i(0)±?(0)ii(t)), 

J Li = (»7L(t)xi(t)»/i(0)ii(0)), (24) 
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FIG. 1: (Color online) Plots of [\g (t)]\ vs. i for TV = 4 and TV = 
8. The parameters for the figure are Tl ~ 2.0, Tr — 2.0, k = 
1.0, k = 0.0, k' = 0.0, 7l = 7_r = 2.5 and A = 0.4. Here 
[[go (£)]] denotes disorder averaged go(t). The average is done over 
100 disorder realisations. Inset shows the plots of Ai(t) and AN(t) 
for TV = 8 for a single disorder configuration. 



where t > 0. 

Now we will calculate all these J's using Eq. (13) and Eq. (5). 
We will present the calculation of Jli explicitly and state the 
results for other J's. Putting the form of x\(t) in the expres- 
sion of Jli in Eq. (24) we get 

/t r-t M M 

dh I dt 2 I dt 3 / dt 4 
- OO J — OO J — OO J — OO 

x K 1 (t 1 ,t 2 ,t 3 ,t 4 ,t), (25) 

Where t 2 , t 3 , t^, t) is given by 

K\ (ti, t 2 , t 3 , *4, t) 

[[Gn(t - h^h) + G 1N (t - t^mih)] 

x[G n {t - t 2 )rj L (t 2 ) + G 1N (t - t 2 )r] R {t 2 )] 
x [Gu(-t 3 )»7i(t3) + G 1N (-t 3 ) m (t 3 )} 
x[G n (-h)r] L (t 4 )+G 1N (-t 4 )ri R (t 4 )}) . (26) 

After taking the average over noises and using their Gaussian 
property, we get 

Ki(t\,t 2 ,t 3 , t 4 , t) 

= A{K^\ti,t 2y t 3 , U, t)S(ti - t 2 )S(t 3 - t 4 ) 



+K[ 2 \t 1 ,t 2 ,t 3 ,h,t)S(t 1 -t 3 )5(t 2 -h) 
+K<f ) (t 1 ,t 2 ,t 3 ,t 4 ,t)5(t 1 -t 4 )5(t 2 ~t 3 )) 



(27) 



where expressions for these K-y are given in Appendix B. 
Putting the expression of K\(t\,t 2 ,t 3 ,t 4 ,t) in Eq. (25) and 
arranging the terms we finally get 



Jli = ±iI[{ilT l A 1 {0)+ 1r T r A n {0)Y 
+ 2{ lL T L A 1 {t)+ lR T R A N {t)} 2 }, 



(28) 



where we have used the definitions of Ai (t) in Eq. (15). Sim- 
ilarly we calculate other J's and their expressions are 

Jl2 = -AtHTl [^Gii(O) { 1L T L A! (0)+7 fi T HJ 4jv(0)}, 

Jls = -4 7 iT L [iG 11 (0){ 7L r L A 1 (0)+ irT r A n (0)} 

+ 2G 11 (t){ lL T L A 1 {t) + lR T R A N (t)}}, 
Ju = ^lLT L W){lLT L A l {t)+ lR T R A N {t)} 




N=128 

N=256 

N=512 

N=1024 

N=2048 

N=4096 

N=8192 



(29) 



Collecting all the expressions for J's from Eqs. (28) and (29) 
in Eq. (24) and subtracting (j'l) 2 we finally obtain 

C AT (t) = 4 7L T L { 7L T L A 1 (0) + lR T R A N (0)}S(t) 

- %li 2 [{lLT L A 1 {t)+ lR T R A N (t)} 

x {T^Mt) + (2T L - T R )j R A N {t)}] , 

= ±lLT L { 1 LT L A 1 (0)+ lR T R A N {Q)}5{t) 

- 9A T (t) (30) 

where 

g AT ( t ) = S lL 2 [{ lL T L A 1 {t)+ lR T R A N {t)} 

x {T LlL A 1 (t) + (2T L -T R ) lR A N (t)}] (31) 

and we have used the identity in Eq. (16). From the above 
expression of g/vr{t) we note that go(t) is always positive. 
Thus we have obtained a closed form expression for the 
non-equilibrium current-current auto-correlation function ex- 
pressed in terms of the Green's function for a disordered har- 
monic chain of length N. The delta function appearing in 
the above equation is purely due to the white nature of the 
noises. More generally one can define the current operator on 
any bond on the harmonic chain. However the detailed form 
of the bond-correlation function is quite different from that of 
the boundary-correlation function. The notable difference that 
we find is the absence of the <5-function peak. We have veri- 
fied that the integral of bond-correlation agrees with the value 
for the boundary-correlation. 



IV. Discussions 

In this section we plot the function gArrif) = 
±lLT L { lL T L A 1 {0)+ lR T R A N {0)}8{t) - C AT . To find the 
functional form of gArif) we need to know the functional 
forms of the functions Ai(t). These functions can be obtained 
by Fourier transforming lu 2 \Gu(ijj)\ 2 as shown in Eq. (19). 
For a general N-particle mass disordered chain it is difficult 
to find analytical expressions for the functions |G^ (w)| 2 . For 
the ordered case G^-(w) can be obtained analytically using 
the tridiagonal nature of the force matrix $ (see for example 
Refn. [12]). However in case of disordered chain, Gf N (to) 
and Grjj (to) can be obtained through transfer matrix approach 



FIG. 2: (Color online) Plots of [[Ajv(t)]] vs. t for different system 
sizes. The parameters for the figure are same as those for Fig. 1. 
A = 0.4. 



in which Gf N (to) and G^(oi) are expressed in terms of a 
product of N random matrices [11]. We numerically evalu- 
ate G^ n (uj) and G^cj) using this transfer matrix approach. 
We observe that at large to > tu d = j^, [\\G^ N (u))\ 2 T\ de- 

,-aNui 2 



cays as e -1 ""^ (a is a positive constant) where m = [[mi]] 
and a 2 = [[(mi — to) 2 ]]. Here [[...]] denotes disorder average. 
This behaviour was proved analytically by Matsuda and Ishi 
[15] and was first observed numerically by Dhar [11]. Another 
observation made by Dhar was that for ui < uJd disordered av- 
erage of G^ n (uj) | 2 is almost identical to that of an ordered 
chain for both the BCs. We make use of this observation in 
this paper. Another observation which we made is that for 
cj > u m the function |Gii(u>)| 2 decays as 1/uj 4 , where u> m is 
the maximum normal mode frequency. This 1 /w 4 behaviour 
can be easily obtained through the transfer matrix approach. 
For small frequencies disorder average of |G^ (u>) \ 2 oscillates 
with u and is again identical to that of ordered chain. 

After integrating Eq. (19) numerically, we obtain A, (t) and 
Gij it) and hence go (t) for different system sizes with differ- 
ent disorder configurations. In Fig. 1 we plot [[<?at(£)]] versus 
t for system sizes N = 4, 8 and 16 with free BC. We observe 
that the correlation functions for two system sizes remains al- 
most identical at short times and starts being different signifi- 
cantly after some time scale. These observations can be made 
by looking at the dominant contributions of w 2 |G^(o;)| 2 in 
the integrand of Eq. (19) for fixed t. At large uj the functions 
\Gf N (u>)\ 2 decays as e~ aN ^ (a is apositive constant)[l 1, 15] 
whereas ^^(a;)! 2 decays as 1/u; 4 . At small frequencies both 
Gi N (uj) and G^ (u>) are oscillating function of us and the fre- 
quency of oscillation increases with system size N. As a re- 
sult A\ (t) is independent of system size N at small times and 
starts depending on N after some time scale, where contri- 
bution from small to becomes important. Whereas, in case 
of A]si(t), only a small range of to contribute in the Fourier 
transform of oj 2 \G^ N (oj)\ 2 (Eq. (19)). For large AT, at small 
times Ai(t) is much larger that Ajv(t) and contributes most 
in go(t), which makes go(t) to be independent of N at small 
times. Inset in Fig. 1 compares A\(t) and Ajv(t) for = 8. 



In the next paragraph we will see that physically interesting 
quantities like current, fluctuations in current in NESS are re- 
lated to the time integral of Cat(£) an d this integral depends 
only on An it), though A\(t) has dominant contribution in 
the correlation function itself. Hence it is more relevant to 
see the behaviour of An it) with system size N. In Fig. 2 we 
plot [A jv it)} for different system sizes. Here we prefer to give 
plots of disordered averaged quantities, since very often we 
are interested in disorder averaged quantities. 

Let Q(t) = J" T dtjj^it) be the heat transfer in duration r 
from left reservoir to the system. Using stationarity property 
of the correlation function it is easy to show that the 2 nd order 
cumulant of Q(t) is related to Cat it) as 
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dtC AT it). 



(32) 



Now integrating the expression of Cat(0 given in Eq. (30) 
from to oo and again using the identity in Eq.(16) we get 

/ dtC AT it) = 2 lLlR T L T R A N iO) 
Jo 

/>00 

+ 8jhUTL-T R ) 2 / dtA 
Jo 



Nit) 



(33) 



In the frequency space the Eq. (33) can be written as an inte- 
gration over io of the transmission coefficient T{(o) defined in 
Eq. (22) and we obtain 



dtC^rit) 



;/,. l lt f ' • 



4tt ,,o 



TlTr 
2vr 



dojT 2 iuj) 
duo T(u>) . (34) 



This expression matches with the expression given in [16] for 
quantum mechanical systems, in the high temperature limit. 
Now if we put Tl = Tr = T in the expression in Eq. (33) 
and use Eq.(17) we get a relation between the current in the 
non-equilibrium steady state and the equilibrium correlation 
function similar to the GK relation derived in [7] 



dtCoit) 



2n 



dio T(w) = T 2 



3 



in - t r 



-, (35) 



where Cq it) is the equilibrium auto-correlation function for 
the open system. The inset of Fig. 3 shows the system size 
dependence of the disorder average of current. 

In general for large system sizes [[j]] and scales 
with N as and N~ a respectively. Using the frequency 
dependence of T(w) = [[T(w)]] an d [[T 2 (u;)]] one can pre- 
dict the value of a and (3 for different BC's. By computing 
[[j]} in NESS, several authors have already studied asymp- 
totic size dependence of [[j]}. Rubin and Greer [17] obtained 
(3 = 1/2 for free BC, which was latter proved rigorously by 
Verheggen[18]. Casher and Lebowitz [13] studied the same 
model and obtained a lower bound for [[j]] > N~ 3 ^ 2 and 
simulations by Rich and vischer [19] confirmed the exponent 



FIG. 3: (Color online) This figure shows the dependence of non- 
equilibrium current fluctuation on system size for free BC. The pa- 
rameters for the figure are same as those for Fig. 1 except Tl = 3.0 
and Tr — 2.0. Inset shows the dependence of non-equilibrium cur- 
rent on system size for free BC. Disorder average is taken over 100 
different disorder realizations. Standard deviation corresponding to 
each point is smaller that the size of the point symbol. 



to be (3 = 3/2. Later Dharfll] obtained j for both the bound- 
ary conditions using Langevin Equation and Green Function 
approach and obtained (3 = 1/2 for free BC and (3 = 3/2 
for fixed BC. Here we follow the same procedure described 

in [11] to find the asymptotic size dependence of [[ — ^"^ c ]] 
from the expression given in Eq. (34). 

We numerically observe that for both the BCs [[T 2 (w)]] is 
much smaller than T(w) for each N. Hence, in determining 
the assymptotic N dependence, dominant contribution comes 
from the integration of T(uj) over uj. To determine a, we use 
the fact (discussed in the first paragraph of this section) that 
for to greater than u>d ~ iV -1 / 2 , T{io) decays exponentially 
as e~ aNu whereas, for u> < LJd, T(w) is almost identical 
to Toioj) of an ordered chain. It can be shown that trans- 
mission coefficient of an ordered chain, denoted by Toiui), 
is independent of u for free BC and goes as w 2 for fixed 
BC. Now putting these forms of 7^(w) and integrating up- 
to io d - N~ 1/2 we get a = 1/2 for free BC and 3/2 for 
fixed BC. We see that the asymptotic size dependence of cur- 
rent fluctuation is same as that of NESS current. We numer- 
ically evaluate the RHS of Eq. (33) for free BC and obtain 
^ ^ a for t —> oo for different system sizes. In Fig. 3 we 
plot [[iS-fc^s]] versus system size N, which shows that the 
fluctuation in current scales with system size as N~ x / 2 , when 
both ends of the chain are free. In the pinned case, since there 
are no low frequency modes, T(w) decays exponentially and 
hence fluctuations in current decays exponentially with N. 



V. Conclusion 

In conclusion we have given an expression for the current- 
current correlation for a one dimensional mass-disordered har- 
monic system in NESS. The correlation function has been ex- 



6 



pressed in terms of the phonon Green's functions which are 
easy to evaluate numerically. We show that the integration of 
equilibrium correlation function gives current satisfying the 
finite size open system Green-Kubo formula whereas the in- 
tegration of non-equilibrium correlation function gives infor- 
mation about current fluctuation in the NESS. Using the non- 
equilibrium correlation function we obtain asymptotic system 
size scaling of the fluctuation in the steady state current. A 
possible application of our results is that they can serve to 
test numerical codes for simulations studying correlations in 
non-equilibrium systems. In this paper we have considered a 
classical one dimensional system with white noise Langevin 
dynamics. It will be straightforward to get an expression for 
correlation function in quantum systems and higher dimen- 
sional systems. 



Appendix A: Proof of Eq.[16] 

Let us first define few quantities: 

G = M^GM? 
f = M-JTM-5 

Using this above definitions Eq. (11) can be written as 

G{t) + fd(t) + *G(t) = (Al) 



We use the above equation to evaluate -^[G T (t')G(t' + t)] 



and get 

dt 



d 7 [G T {t')6(t' +t)] = -2G T {t')tG{t' + t) 
+ ^lG~T(t>)<bG(t' + t)] 



Now integrating both side of the above equation over t' = 
to t' = oo we get 



G{t) = 2 / dt'G 1 (t')TG(t' +t). 



(A2) 



To the above equation we have used the following: G(0) = 
M-\ G(0) = 0, G(t) ->■ as t ->■ oo. Now we know that 
r« = {%5 a + ^5 iN )5ir Taking (ll)*' 1 element on the 



both side of the matrix equation (A2) we get 
Gu(t) 



poo 

/ df[7 
Jo 



L Gu(t')Gu(t' + t) 



(A3) 



+ j R G 1N (t')G 1N (t' + t)} 
= jLAxit) + j R A N (t) 

Appendix B: Expressions of K\ 's 

[l 2 L TlGn(t - h)G n (t - t^G^-t^G^-U) 

+7i?7 1 R,GiAr(t — t\)GiN{t - *2)Giat(— i3)GlJv(— *4) 

+l L T LlR T R {G 1N {t - h)G 1N (t - t 2 )G 11 (-t 3 )G 11 (-U) 
+Gn(t - ti)Gn(t - t 2 )G 1N (-t 3 )G 1N (-U)}] , 



K[ 2 \t 1 ,t 2 ,t 3 ,U,t) = 

[llTlG n (t - h)G n (t - i 2 )G 11 (-i 3 )G 11 (-i 4 ) 

+7fl^RG'iAr(i — t\)G\N{t - t 2 )GiN{—t 3 )GiN{—ti) 

+jLT L j R T R {G 1N (t - t!)Gu(i - t 2 )G 1N (-t 3 )G 11 (-U) 
+G n (i - ti)Gi N (t - t 2 )Gu(-t 3 )G 1N (-U)}] 



and 



K{ 3) (h,t 2 ,t 3 ,U,t) = 

[llTlGuit - ti)Gu(t - t2)Gu(-t3)Gu(-U) 
+ lR T 2 R G 1N {t - h)G 1N (t - t 2 )G 1N (-t 3 )G 1N (-U) 
+lLT LlR T R {G 11 {t - h)G 1N (t - t 2 )G 1N (-t 3 )Gn(-t4) 
+G 1N {t - ti)Gn(t - t2)Gn(-t3)Gij V (-t 4 )}] 



I thank Dr. Abhishek Dhar for useful suggestions and Jayaku- 
mar A for helpful discussions. 
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